---
title: "R Markdown for OLS with GETS procedure & Correspondence Analysis"
author: "Ulrich Fritsche & Gabi Waldhof"
date: "2022-11-24"
output:
  html_document: default
  word_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Idea/ Structure

1. Install necessary packages
2. Activate necessary libraries
3. Load the data
4. Estimate a linear probability model/ tests
5. Apply HEC standard errors 
6. Test for Multicollinearity with Variance Influence Factor
7. Use GETS procedure
8. Estimate model based on GETS results
9. Stargazer output
10. Correspondence Analyses of Narratives about GE

## Install packages

```{r}
#install.packages("readxl")
#install.packages("parallel")
#install.packages("zoo")
#install.packages("lmtest")
#install.packages("estimatr")
#install.packages("car")
#install.packages("stargazer")
#install.packages("gmodels")
```

## Activate libraries

```{r}
library(readxl)
library(gets)
library(summarytools)
library(lmtest)
library(estimatr)
library(sandwich)
library(gets)
library(car)
library(stargazer)
library(gmodels)
library(MASS)
```

## Load data

```{r}
DF.OLS1 <- read_excel("data_Narratives about GE_cleaned and transformed.xlsx")
DF.OLS1["Degree"][DF.OLS1["Degree"] == 7] <- NA
```

## Estimate OLS and make tests

```{r}
MOD.OLS1 <- lm(GeAttitude ~ Gender + Age + Region + Income + Degree + Interest + Importance + Meaning + Voting + Sentiment1 + Sentiment2 + Sentiment3 + Sentiment4 + Sentiment5 + Emotionality + EmoContempt + EmoAnger + EmoFear + EmoDisgust + EmoJoy + EmoSurprise + EmoSadness + RegStrict + RegFieldUni + RegFieldComm + RegLabUni + RegLabComm + RegImport + RegExport + RegChoice + RegChoiceLabel + RegNoChoice + MERGED_MFT_NARRATIVE1, data = DF.OLS1)
summary(MOD.OLS1)
coeftest(MOD.OLS1, vcov = vcovHC(MOD.OLS1, type = "HC3"))
OLS1.res = resid(MOD.OLS1) # save resids
plot(OLS1.res) # plot resids
hist(OLS1.res) # plot distribution
summarytools::descr(OLS1.res) # show descriptive stats
lmtest::bptest(MOD.OLS1) # Breusch-Pagan test
```

## Apply HEC, compare with stargazer

```{r}
# Adjust standard errors
cov1         <- vcovHC(MOD.OLS1, type = "HC3")
robust_se1    <- sqrt(diag(cov1))

# Stargazer output (with and without RSE)
stargazer(MOD.OLS1, MOD.OLS1, type = "text",
          se = list(NULL, robust_se1))
```


## Test for Multicollinearity with Variance Influence Factor

```{r}
vif(MOD.OLS1)
```

## Apply general-to-specific-procedure (gets)

```{r}
#ISAT <- isat(MOD.OLS1, vcov.type = "white")
#GETS.ISAT.MOD.OLS1 <- gets.isat(ISAT)
#MOD.GETS <- as.lm(GETS.ISAT.MOD.OLS1)
#summary(MOD.GETS)
```



```{r echo=T, results='hide'}
GETS.MOD.OLS1 <- gets.lm(MOD.OLS1)
```


## Estimate the model with the results general-to-specific-procedure (gets)

```{r}
formula.ols2 <- GeAttitude ~ Gender + Age + Region + Income + Degree + Interest + Importance + Meaning + Voting + Sentiment1 + Sentiment2 + Sentiment3 + Sentiment4 + Sentiment5 + Emotionality + EmoContempt + EmoAnger + EmoFear + EmoDisgust + EmoJoy + EmoSurprise + EmoSadness + RegStrict + RegFieldUni + RegChoiceLabel + MERGED_MFT_NARRATIVE1 
MOD.OLS2 <- lm(formula.ols2, data = DF.OLS1)
summary(MOD.OLS2)
coeftest(MOD.OLS2, vcov = vcovHC(MOD.OLS1, type = "HC3"))

# Adjust standard errors
cov2         <- vcovHC(MOD.OLS2, type = "HC3")
robust_se2    <- sqrt(diag(cov2))

# Stargazer output (with and without RSE)
stargazer(MOD.OLS2, MOD.OLS2, type = "text",
          se = list(NULL, robust_se2))
```

## Stargazer Output

```{r}
stargazer(MOD.OLS1, MOD.OLS2, type = "text", 
          se = list(robust_se1, robust_se2))
```


## Correspondence Analyses of Narratives about GE

```{r}
library(readxl)
#install.packages("gmodels")
library(gmodels)
```

#### 1st & 2nd Narrative

```{r}
Correspondence_Analysis <- read_excel("correspondence_1_2.xlsx")
CrossTable(Correspondence_Analysis$MERGED_MFT_NARRATIVE1, Correspondence_Analysis$MERGED_MFT_NARRATIVE2, chisq = TRUE, expected = TRUE, sresid = TRUE, format = "SPSS")
```

#### 1st & 3rd Narrative

```{r}
correspondence_1_3 <- read_excel("correspondence_1_3.xlsx")
CrossTable(correspondence_1_3$MERGED_MFT_NARRATIVE1, correspondence_1_3$MERGED_MFT_NARRATIVE3, chisq = TRUE, expected = TRUE, sresid = TRUE, format = "SPSS")
```

### 2nd & 3rd Narrative

```{r}
correspondence_2_3 <- read_excel("correspondence_2_3.xlsx")
CrossTable(correspondence_2_3$MERGED_MFT_NARRATIVE2, correspondence_1_3$MERGED_MFT_NARRATIVE3, chisq = TRUE, expected = TRUE, sresid = TRUE, format = "SPSS")
```